High-frequency dynamical response of Abrikosov vortex lattice in flux-flow region 
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The dynamical response of the Abrikosov vortex lattice in the presence of an oscillating driving 
field is calculated by constructing an analytical solution of the time-dependent Ginzburg-Landau 
equation. The solution is steady-state, and work done by the input signal is dissipated through 
vortex cores, mainly by scattering with phonons. The response is nonlinear in the input signal, and 
is verified for consistency within the theory. The existence of well-defined parameters to control 
nonlinear effects is important for any practical application in electronics, and a normalised distance 
from the normal-superconducting phase-transition boundary is found to be such a parameter to 
which the response is sensitive. Favourable comparison with NbN experimental data in the optical 
region is made, where the effect is in the linear regime. Predictions are put forward regarding the 
suppression of heating and also the lattice configuration at high frequency. 
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I. INTRODUCTION 

Superconducting films are candidate substances for the 
improvement of electronics technology in a myriad of ap- 
plications. While the low resistance is very attractive in 
this regard, it has proved difficult to control the nonlin- 
ear behaviour of such materials in response to electro- 
magnetic field 1 . When a magnetic field is strong enough 
to penetrate into a superconductor in the form of quan- 
tised magnetic flux tubes, the vortex state obtains as a 
mixed state of superconducting phase punctuated by the 
vortices themselves. Vortices are surrounded by a su- 
percurrent and can be forced into motion by the current 
resulting from an applied field. 

As a topological defect, a vortex is not only stable 
under perturbations^ but cannot decay. The collec- 
tion of vortices in a type-II superconductor forms what 
is called vortex matter, and it is this which determines 
the physical properties of the system rather then the un- 
derlying material properties, in particular driving phase 
transitions[2.4]. In the mixed state, a superconductor is 
not perfect; it exhibits neither perfect diamagnetism nor 
zero electrical resistance. The transport current J gen- 
erates a Lorentz force F = $oJ x h on the vortex and 
forces it into motion, dissipating energy. 

In reaching thermal equilibrium, energy is transferred 
via interactions between phonons and quasiparticle exci- 
tations. Small-scale imperfections such as defects scatter 
the quasiparticles, affecting their dynamics. In dirty su- 
perconductors, impurities are plentiful and vortices ex- 
perience a large friction. This implies a fast momentum- 
relaxation process. In contrast is the clean limit, where 
impurities are rare and no such relaxation process is avail- 
able. It is in this situation of slow relaxation that the Hall 
effect appears. 

Generally, the H-T phase diagram^ of the vortex mat- 
ter has two phases. In the pinned phase vortices are 



trapped by an attractive potential due to the presence of 
large-scale defects, thus resistivity vanishes. This phase 
contains what are known as glass states. There is then 
the unpinned phase in which vortices can move when 
forced and so a finite resistivity appears. This phase 
is also known as the flux-flow region and can be of two 
types. One type is a liquid state where vortices can move 
independently; the other type is a solid state in which 
vortices form a periodic Abrikosov lattice* resulting from 
their long-range interacton. One model for the transition 
between the pinned and unpinned phases appears in [f|. 

In the unpinned phase, the system is driven from equi- 
librium and experiences a relaxation process. There 
are several ways to describe such a system. A micro- 
scopic description^ invoking interactions between a vor- 
tex and quasi-particle excitations at the vortex core pro- 
vides a good understanding of friction and sports good 
agreement with experiments in the sparse-vortex region 
H <C H c2 . There is also a macroscopic description, the 
London approach, where vortices are treated either as in- 
teracting point-like particles or an elastic manifold sub- 
ject to a pinning potential, driving force and friction 8,9 . 
In the small-field region, vortices behave as an array of 
elastic strings. 

In the dense- vortex region H 3> H c \, where the mag- 
netic field is nearly homogeneous due to overlap between 
vortices, Ginzburg-Landau (GL) theory, which describes 
the system as a field, provides a more reasonable model. 
In dynamical cases, time-dependent GL (TDGL) theory 
is appropriat ed 10 ' 1 1 ; in GL-type models, additional sim- 
plification can come from the lowest Landau level (LLL) 
approximation which has proven to be successful in the 
vicinity of the superconducting-normal (S-N) phase tran- 
sition line H ~ H c i- This has been pursued in the static 
case^ (without driving force) and in the dynamic case 
with a time-independent transport current^. It may be 
noted that in the glass state, zero resistance within the 
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LLL approximation cannot be attained 13 . 

Based on TDGL theory, we will study the dynamical 
response of a dense vortex lattice forced into motion by 
an alternating current induced by an external electro- 
magnetic field. Vortices are considered which are free 
from being pinned and thermally excited, which in addi- 
tion to thermal noise would produce entanglement and 
bending. We assume the vortices can transfer work done 
by an external field to a heat bath. Experimentally, a 
low-temperature superconductor far away from the clean 
limit is the best candidate for attaining these conditions. 
We do not consider thermal fluctuation effects specific to 
high-temperature superconductors. In a dissipative sys- 
tem driven by a single-harmonic electric field E cos £It, 
long after its saturation time we can expect the system 
to have settled into steady-state behaviour, where the 
vortices are vibrating periodically with some phase. 

The TDGL model in the presence of external electro- 
magnetic field is analysed and solved in section [Til The 
dynamical S-N phase transition surface T C (H, E, f2) is lo- 
cated in {T, H, E, 57}-space. This surface coincides with 
the mean- field upper-critical field H C 2(T) in the absence 
of the applied field, and with the phase-transition surface 
in the presence of the constant driving field considered 
by Hu and Thompson 1 ^. We will provide an analytical 
formalism for perturbative expansion in the distance to 
T c , valid in the flux-flow region. The response of vor- 
tex matter forced into motion by the transport current is 
studied in section IIIII The current-density distribution 
and the motion of vortices are treated in section IHI Al 
In analysing the vortex lattice configuration in section 
IIIIB1 a method is utilised whereby the heat-generation 
rate is maximised. Next are discussed power dissipa- 
tion, generation of higher harmonics, and the Hall effect. 
An experimental comparison is made in section IIVI with 
Far-Infrared (FIR) measurement on NbN. Finally, some 
conclusions are made in section [V] 



II. FLUX-FLOW SOLUTION 

Let us consider a dense vortex system prepared by ex- 
posing a type-II superconducting material to a constant 
external magnetic field H = (0,0, —H) with magnitude 
H c2 > H ^> H c ±. We also select the c axis of the su- 
perconductor to be in the z direction. Let the supercon- 
ductor carry an alternating electric current along the y 
direction, generated by an electric field E(t) = Ecosilr 
as shown in Fig. [TJ Such a system when disturbed from 
its equilibrium state will undergo a relaxation process. 
For our system, the TDGL equation 10,14 is a useful ex- 
tension of the equilibrium GL theory. 

In the dense-vortex region of the H-T phase diagram, 
vortices overlap and a homogeneous magnetic field ob- 
tains. Describing the response of such a system by a field, 
the order parameter '5 in the GL approach, is more suit- 
able than describing vortices as particle-like flux tubes, 
as is done in the London approach^. 



A. Time-dependent Ginzburg-Landau model 

A strongly type-II superconductor is characterised by 
its large penetration depth A and small coherence length 
£, k = A/^3>1. The difference between induced mag- 
netic field and external magnetic field is H B = — 47rM. 
In the vicinity of the phase-transition line H C 2(T) vor- 
tices overlap significantly, and H ps B making M small. 
In this case, the magnetic field may be treated as homo- 
geneous within the sample. We will have in mind an 
experimental arrangement using a planar sample very 
thin compared with its lateral dimensions. Since the 
characteristic length for inhomogeneity of electric field 10 
£g = 47rA 2 <7 rl /7 is then typically large compared with 
sample thickness, this implies that the electric field may 
also be treated as homogeneous throughou t 10 ^ 2 , elimi- 
nating the need to consider Maxwell's equations explic- 
itly. 

In equilibrium, the Gibbs free energy of the system is 
given by^ 

F i*i = / dr {^ |M|2 + £: |a =* 12 

-a(7j-T)|*| a + ||*|4 (1) 

where T° is the critical temperature at zero field. Co- 
variant derivatives employed here preserve local gauge 
symmetry and are two-dimensional; D T = d T + *%■ A T 
and D = V< 2 ) - i^-A. 

he 

Governing the dynamics of the field \P is the TDGL 
equation 

h 2 1 SF 

— —D T * = - — . 2 

2m ab 5^* 

This determines the characteristic relaxation time of the 
order parameter. Microscopic derivation of TDGL can 
be found in 0[l5[ in which the values of a, /3 and 7 
are studied. In the macroscopic case, these are viewed 
simply as parameters of the model. At microscopic scale, 
disorder is accounted for by 7, the inverse of the diffusion 
constant; the relation of 7 to normal-state conductivity 
is discussed in appendix [A"l 

In standard fashion, E = —VA T — \d T A while B = 
V x A. Our set of equations is completed 2 by including 
Ampere's law, writing for the total current density 

J o = (t-) V x V x A = ^nE-r-Jb- (3) 

471" 

As we shortly make a rescaling of quantities, we have 
written subscripts here for clarity. The first term is the 
normal-state conductivity. The second term can be writ- 
ten using a Maxwell-type equation relating the vector- 
potential with the supercurrent, 

he* 

Jo = -i [**D* - *(D*)*] . (4) 

2m a b 
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This is a gauge-invariant model; we fix the gauge 
by considering the explicit vector potential A(r) = 
(By, 0, 0) and A T (r, r) = yE cos fir, corresponding to an 
alternating transport current. Each vortex lattice cell 
has exactly one fluxon. We do not assume the electric 
field and the motion of vortices are in any particular di- 
rection relative to the vortex lattice, by way of rendering 
visible any anisotropy. 

For convenience, we define some rescaled quantities. 
The rescaled temperature and magnetic field are t — 
T/T° and b = H/H° 2 . H° 2 denotes the mean field upper- 
critical field, extrapolated from the T° region down to 
zero temperature. 

In the a-b plane of the crystal we make use of magnetic 
length We define Q = \ 2 jb where £ 2 = h 2 /2m ab aT°. 
The scale on the c-axis is £ C /Vt> with £^ = h 2 /2m c aT® . 
The coordinate anisotropy in z is absorbed into this 
choice of normalisation, as can be seen in ([6]). The or- 
der parameter ^ is scaled by y/2baT®/ f3. The time scale 
is normalised as r s = 7^/2. Therefore, frequency is 
u) = Or, . Note that uj is then inversely proportional to b. 
The amplitude of the external electric field is normalised 
with E = 2h/e*£ 3 j so that e = E/E , 

After our rescaling the TDGL equation takes the sim- 
ple form 
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L^> - — (1 - t)^J> + * = 



where the operator L is defined as 



L = D T 



1 



1 



(5) 



(6) 



crystal structure determines linear combinations of these 
basis elements with respect to k x ; the resulting (p n func- 
tions are then useful for expansion purposes below. The 
inner product is (<Pm,<Pn) — i'fimfn) where the brackets 
(• • • ) denote an integral 2 over space and time. To define 
averages over only time or space alone, we write (• • -) T 
or (■ ■ • ) r respectively. 
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FIG. 1: Possible vortex lattice configurations: (a) Typical 
large angle, plotted for 9 a = 60°; (b) Typical small angle, 
plotted for 9 a = 30 . In the static case, both 9a — 60 and 
9a = 30 are two particular angles will correspond to the 
same energy. The applied constant magnetic field H is along 
the — z direction and time-dependent oscillating electric field 
E(t) is in the y direction. 9a is the the apex angle of the two 
defining lattice vectors. The two vectors for the rhombic vor- 
tex lattice are ai = (^/2ty/p, 0) and &2 = (\/2nJp/2, ^2np) 
where p — ^ tan 8 a- The motion of vibrating vortices, indi- 
cated by the red arrow, is disccused in section [ill AI 



With our specified vector potential, covariant derivatives 
are D T = d T + ivy cos lot, D x = d x — iy and D y — d y . We 
define v = eb~^l 2 for convenience. The TDGL equation 
is invariant under translation in z, thus the dependence 
of the solution in the z direction can be decoupled. L is 
not hcrmitian; 



D T 



1 



1 



(7) 



where the conjugation is with respect to the usual inner 
product, defined below. 

We will make extensive use of the Eigenfunctions of L 
and iJ in what follows. The Eigenvalue equation 



L<p„ 



(8) 



defines the set of Eigenfunctions of L appropriate for our 
analysis; this can be seen in appendix [Bl The convention 
is that e„ = £„< when and only when n — n'. Taking 
corresponding 1 Eigenfunctions of to be (pn,k w , the or- 
thonormality ((p n ,k m , ¥>n' ,k'J = S nn >S(k x -k' x ) may be cho- 
sen, so long as (<fi n ,k X ) fn,k x ) 7^ 0. Shown in appendix IB1 



B. Solution of TDGL equation 

States of the system can be parametrised by (t, b, e, lu). 
By changing temperature t, a system with some fixed 
(6, e, ui) may experience a normal-superconducting phase 
transition as temperature passes below a critical value 
t c (b,e,uj). Such a point of transition is also known as a 
bifurcation point. 

The material is said to be in the normal phase when 
* vanishes everywhere; otherwise the superconducting 
phase obtains, with ^> describing the vortex matter. Be- 
cause of the vortices, the resistivity in the superconduct- 
ing phase need not be zero. The S-N phase-transition 
boundary t c (b, e,w) separates the two phases. To study 
the condensate, we will use a bifurcation expansion to 
solve (0. We expand * in powers of distance from the 
phase transition boundary t c . 



1 One 'follows the sign' in front of the D T in L and switches it 2 Shortly we will be dealing with a periodic system, and we will 
in the resulting ip to get the 'corresponding' Eigenfunction (p for normalise such integrations by the unit cell volume and the pe- 

. riod in time. 
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1. Dynamical phase-transition surface 

As in the static case, we can locate the dynamical 
phase-transition boundary by means of the linearised 
TDGL equation 2 *!!. This is because the order param- 
eter vanishes at the phase transition, and we do not need 
to consider the nonlinear term. The linearised TDGL 
equation is written 



L*-^(l-*)* = 0. 



(9) 



Of the Eigenvalues of L, only the smallest one Eq, cor- 
responding to the highest superconducting temperature 
t c , has physical meaning. The S-N phase transition oc- 
curs when the trajectory in parameter space intersects 
with the surface 

eo-^(l-<) = (10) 
where the lowest Eigenvalue is calculated in appendix |B| 

Utilising a 6-independent frequency v = bu and ampli- 
tude e of input signal, we write 

e 2 /2 

'-t-^wh- (12) 

In the absence of external driving field, e = 0, the phase- 
transition surface coincides with the well-known static- 
phase transition line 1 — t — 6 = in the mean-field ap- 
proach. With time-independent electric field at v = 0, 
where the vortex lattice is driven by a fixed direction of 
current flow, the dynamical phase-transition surface co- 
incides with that proposed in (Toj . but with a factor of 
1/2. This amplitude difference is familiar from elemen- 
tary comparisons of DC and AC circuits. 

In the above equation, we can see that in the static 
case e = 0, the superconducting region is 1 — t — b > 0. 
In addition when e ^ 0, the superconducting region in 
the b-t plane is smaller than the corresponding region in 
the static case, as can be seen in Fig. Ufa). 

Finally, increasing frequency will increase the size of 
the superconducting region, as in Fig. [2jb) ; in the high- 
frequency limit, the area will reach its maximum, which 
is the superconducing area from the static case. As with 
any damped system, response is diminished at higher fre- 
quencies. 

The superconducting state does not survive at small 
magnetic field; for example at e = 0.2 in Fig. [2fa), the 
material is in the normal state over most of the H-T 
phase diagram. Later in this paper we will consider in- 
terpretation of this phenomenon. In particular, when 
discussing energy dissipation in section UlI C[ we will see 
that the main contribution to the dissipation is via the 
centre of the vortex core. At small magnetic field, since 
there are fewer cores to dissipate the work done by the 
electric field, the superconducting state is destroyed and 
the order parameter vanishes. 





(a) 
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FIG. 2: Dynamical superconducting-normal phase transition: 
(a) Critical temperature t c as a function of b for various e at 
v — 0.1 and (b) t c as a function of e for various v at 6 = 0.1. 
The straight line in (a) is the e = curve and corresponds to 
the mean-field phase-transition line 1 — t — b— 0. States above 
each line are normal phase while the region below each line 
is superconducting, e suppresses the superconducting phase 
as shown in (a), while v removes this suppression effect, as 
shown in (b). 



2. Perturbative expansion 

That the vortex matter dominates the physical proper- 
ties of the system is especially pronounced in the pinning- 
frcc flux-flow region. Here we solve (JSJ by a bifurcation 
expansio n 12 ' 16 . Since the amplitude of the solution grows 
when the system departs from the phase transition sur- 
face where *S> = 0, we can define a distance from this 
surface as 



2b 



(l-t)-e , 



and expand $ in e. The TDGL in terms of e is 



L* - + *|*| 







(13) 



(14) 



where L = L — eo is the operator L shifted by its smallest 
Eigenvalue. ^ is then written 



(15) 



i=0 



and it is convenient to expand $W j n terms of our Eigen- 
functions of L 



$(*) = 



n=0 



(16) 



In principle, all coefficients c n in (|15|) can be obtained 
by using the orthogonal properties of the basis, which are 
explained in appendix [B] Inserting ^ from equation f| 15[) 
into TDGL equation (T14")) . and collecting terms with the 
same order of e, we find that for i = 



L$<°> = 



and for i = 1 



L$W -$(°) + $(°)|$(°)| 2 =0. 



(17) 



(18) 



5 



For i = 2 

L$( 2 ) - $(!) + c (°) 2 (2$«|$ | 2 + $«*$ 2 ) = o (19) 
and so on. 

Observing (TTT)) . the solution for the equation is 

(20) 



where ipo is a particular linear combination of all Eigen- 
functions with the smallest Eigenvalue. 

The coefficient of e 1 / 2 can be obtained by calculating 
the inner product of tpo with (Ql 



1 



(21) 



In the same way, the coefficient of the next order e 3 / 2 , 
can be obtained by finding the inner product of tpo with 
the i — 2 equation (fT9|) . 



^ = (2c«^o^»ko| 2 ) + c«*^o,^». 

^Po — 



71=1 



(22) 

The inner product of </5 m on (fT5)l gives the coefficient for 
m > 



.(i) - 



/WA) 



and 



The solution of TDGL is then 



(23) 
(24) 



1/2 10 



e^J2 ci n ) ^ + °^ /2 )- ( 25 ) 



In this paper we will restrict our discussion to the re- 
gion near t c where the next-order correction can be dis- 
regarded; 



(26) 



We would like to emphasise that our discussion at 
this order is valid in the vicinity of the phase-transition 
boundary and in particular for a superconducting system 
without vortex pinning. In such a system, vortices move 
in a viscous way, resulting in flux-flow resistivity; no di- 
vergence of conductivity is expected. Our results based 
on (|26|) were calculated at e 1 / 2 order, where only the low- 
est eigenvalue n = of the TDGL operator L makes an 
appearance. 

The next-order correction is at order e 3 / 2 , and there is 
now a contribution from higher Landau levels. From the 
symmetry argument in [l6lll7j . as long as the hexagonal 
lattice remains the stable configuration for the system, 
the next-order contribution comes from the sixth Landau 
level with a factor (eg — ^o) 1 - Even in the putative 
case of a lattice deformed slightly away from a hexagonal 
configuration, the next contributing term is n = 2, since 
in our system the lattice will remain rhombic. 



C. Vortex- lattice solution 

The vortex lattice has been experimentally observed 
since the 1960s and its long-range correlations have been 
clearly observed^ with dislocation fraction of the order 
10 -5 . Remarkably, the same techniques can be used to 
study the structure and orientation of moving vortex lat- 
tice with steady current^, and with alternating current 
in the small-frequency regime 3 -. In this subsection, we 
will discuss the configuration of the vortex lattice in the 
presence of alternating transport current in the long-time 
limit. 

In the dynamical case, the presence of an electric field 
breaks the rotational symmetry of an effectively isotropic 
system to the discrete symmetry y — > —y. In contrast, 
a rhombic lattice preserves at least a symmetry of this 
kind along two axes, and the special case of a hexagonal 
lattice preserves sixfold symmetry. 

The area of a vortex cell is determined by the quantised 
flux in the vortex, which is 2-7T in terms of our rescaled 
variables. As shown in Fig. [TJ we choose a unit cell C 
defined by two elementary vectors ai and a2. We will 
first construct a solution for an arbitrary rhombic lattice 
parameterised by an apex angle 9a- 

Consideration of translational symmetry in the x di- 
rection leads to the discrete parameter ki = 2irl/ai = 
y/2wpl. In appendix 151 we show that in the long-time 
limit the lowest-eigenvalue steady-state Eigenfunctions of 
L must therefore combine to form 

oo 

ipo = \ / 2^ J2 e^' (i - 1) e*( a: - , ' sin " T /")w fel (y,T). (27) 



Here ipo is normalised as 

The function u is given by 

u kl (y,T) = c(r)e 2 L y v >i , 

with 

■ [sin 2 0+cos2(ut-0) + ^j sin2(wT-0)l 



c(r) = e 



(28) 
(29) 



In analogy with a forced vibrating system in mechanics, 
a phase 6 = tan _1 aj and a reduced velocity v = vcosd 
have been introduced for convenience in (|29|) . The 
zero electric-field limit, large-frequecy limit and zero fre- 
quency limit are consistent with previous studies con- 
cluded in appendix [Bl 

In our approximation, the /3q in (|26[) is a time- 
independent quantity from (|24[) and 

q = A/£ d(0JT) (e 4( 1 + cos2 " 1 - + ("- 1 /")sin2 W r) 

27T J \ 

oo 



/ !■■ COS CUT 



q— — oo 



(30) 
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where = v/(l + uj 2 ). 

In the small signal limit v — > 0, ph reduces to the 
Abrikosov constant. The Abrikosov constant with cither 
6 a = 30 ° or 60 ° minimises the GL free energy ([TJ in the 
static stated. To be more explicit, /3q can be expanded 
in terms of the amplitude of input signal. In powers of v, 

| = 1 + ^<f-f> + °<« 3 > < 31 > 

and we find it convenient to write in terms of w w . The 
first term in po is the Abrikosov constant 

oo 

p A = ^ £ (-)W e -*(*5+*S>. (32) 

p.q— — oo 

For hexagonal lattices (3a ~ 1.159, whereas for a square 
lattice Pa ~ 1.18. The next term in /3 is v 2 , with a 
coefficient 

y OO 

h = *T £ (-r(^ + ^)e-i^+ fc ?). (33) 

— oo 

We see that at high frequency, the correction in higher 
order terms of v u can be disregarded. 

III. RESPONSE 

In this section we discuss the current distribution and 
motion of vortices, energy transformation of the work 
done on the system into heat, nonlinear response and 
finally the Hall effect. 



A. Motion of vortices 

In addition to the conventional conductivity at- 
tributable to the normal state, there is an overwhelm- 
ing contribution due to the superconducting condensate 
in the flux-flow regime, tempered only by the dissipa- 
tive properties of the vortex matter. In this section we 
will examine the supercurrent density to investigate the 
motion of the vortex lattice. We consider a hexagonal 
lattice in a fully dissipative system; the non-dissipative 
part known as the Hall effect will be discussed in section 

HheI 

The supercurrent density ^7(r, r) is obtained by sub- 
stitution of the solution (j2"6")l into 

j x =j Yu (*4^-»W*^). (34) 

p,g— — oo ^ ' 

and 

Jv = j- £ ( l{kp 2 kq) -g«*(^-g)V,(35) 

' ® p,q— — oo ^ ' 



where 

„ — „-if (p 2 -1 2 -P+l)p i ( k q- k p){ x -T3 sinwr) * 

yp.q — <- e a k p a k q i 

and u is given in (f2"5|). Observing Fig. [31 we conceptually 
split the current into two components. One part is the 
circulating current surrounding the moving vortex core 
as in the static case; we refer to this component as the 
diamagnetic current. The other part which we term the 
transport current is the component which forces vortices 
into motion. 3 




FIG. 3: Current flow at v = 1, w = 1: (a) t = 0; vortex cores 
move to the right, (b) t — 4tt/5; vortex cores move to the 
left. Vortices are drawn back and forth as the direction of 
the transport current density alternates. The magnitude of 
the current density has maximal regions which tend to cir- 
cumscribe the cores; the maxima in these regions move in the 
plane and their manner of motion can be described as lead- 
ing the motion of the vortices by a small phase. The average 
current in a unit cell leads the motion of the vortex in time 
by a phase of n/2. 

The diamagnetic current may be excised from our con- 
sideration by integrating the current density over the unit 
cell C; that is, we consider (J) T . We have (J X ) T = and 

(J) r = -±6cos(wT-0)e-£ sin2 (" T - e ) + ^. (36) 
Po 

With our conventions, the transport current is along the 
y-direction. Considering the Lorentz force between the 
magnetic flux in the vortices and the transport current, 
we expect the force on the vortex lattice to be perpen- 
dicular. 

We identify the locations of vortex cores to be where 
1^1 2 = 0. The velocity of the vortex cores turns out to 
be 

v c {t) = vcos(ujt — 9) (37) 

along the cc-direction. Note that vortex lattice moves 
coherently. The vortex motion the electric field with a 



3 Thinking of the system being embedded in three-dimensional 
space, the circular and transport currents are essentially the curl 
and gradient components of the current. 
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phase 9 which increases with frequency and reaches 7r/2 
asymptotically. The maximum velocity of vortex motion 
v decreases with increasing frequency. 

In Fig. [3] we show the current distribution and the re- 
sultant oscillation of vortices. As anticipated, the trans- 
port current and the motion of vortices (|37[) are perpen- 
dicular as the vortices follow the input signal. The cur- 
rent density diminishes near the core; it is small there 
compared to its average value. 

In steady-state motion, since the vortices move coher- 
ently in our approximation, the interaction force between 
vortices is balanced as in the static case. Since the sys- 
tem is entirely dissipative, the motion that the vortices 
collectively undergo is viscous flow. The vortex lattice 
responds to the Lorentz driving force as a damped oscil- 
lator, and this is the origin of the frequency-dependent 
response. 

B. Configuration of moving vortex lattice 

In static case the system is described by the GL equa- 
tion. Solving this equation, which is ([2]) but with zero 
on the left-hand side, will select some lattice configura- 
tion. The global minima of the free energy correspond 
to a hexagonal lattice, while there may be other con- 
figurations producing local minima. In the static case 
the lattice configuration can be determined in practice 
by building an Ansatz from the linearized GL solution 1 ® 
and then using a variational procedure to minimise the 
full free energy. 

In the dynamic case, there is no free energy to min- 
imise; we must embrace another method of making a 
physical prediction regarding the vortex lattice configu- 
ration. Let us follow [12| and take as the preferred struc- 
ture the one with highest heat-generation rate. Though 
we have at present no precise derivation, our physical 
justification of this prescription is that the system driven 
out of equilibrium can reach steady-state and stay in con- 
densate only if the system can efficiently dissipate the 
work done by the driving force. Therefore, whatever the 
cause, the lattice structure most conducive to the main- 
tenance of the superconducting state will correspond to 
the maximal heat generation rate. 

The heat-generating rate^ is 

(Q) r = 2( \D T iJ>\ 2 ) r (38) 

~2 

_ J_^_ p \ cos 2 0-f- sin2(wr-0) 

A) 4 

v 2 

{ cos2(o;t - 6) + 1 + — [cos4(o;r - 9) + l]}. 
8 

fin is given explictly in (|30[) and is the only parameter 
involving the apex angle 9a of the moving vortex lattice. 
Here /3o plays the same role as the Abrikosov constant ft a 
in the static case. Corresponding to maximising the heat- 
generating rate, the preferred structure can be obtained 
by simply minimising /3q with respect to 9a- This shows 



from the current viewpoint of maximal heat-generation 
rate that vortices are again expected to move coherently. 

In ([31})) or (f3"T|) it is seen that the moving lattice is 
distorted by the external electric field but this influence 
subsides at high frequency. Numerical solution for min- 
imising 9a shows that while near the high-frequency limit 
there remain two local minima for /3o with respect to 9a, 
the solution near 9a = 60 is favoured slightly over that 
at 30 as the global minimum. This is as presented in 
Fig. [TJa). The two minima tend to approach each other 
slightly as the frequency begins to decrease further. In 
an experimental setting, this provides an avenue for test- 
ing the empirical validity of the maximal heat generation 
prescription, in particular in terms of the direction of 
lattice movement^. 

We put forth the physical interpretation that at high 
frequency the friction force becomes less important, and 
the distortion is lessened. Since interactions dominate 
the lattice structure the system at high frequency will 
have many similarities with the static case. 



C. Energy dissipation in superconducting state 

Energy supplied by the applied alternating current is 
absorbed and dissipated by the vortex matter, and the 
heat generation does not necessarily occur when and 
where the energy is first supplied. In Fig. 01 we show 
an example of this transportation of energy by the con- 
densate. On the left is shown a contour plot of the 
work (P) T — (J ■ v) r done by the input signal; points 
along a given contour are of equal power absorption. On 
the right of Fig. [4] is shown the heat-generating ratei^, 
(Q)t = 2(|-DtV| 2 )t- The periodic maximal regions are 
near the vortex cores in both patterns. 




FIG. 4: Work contours of superconducting component at 
v — 1 and uj = 1: (a) Work (P) T and (b) Heat generating 
rate (Q) T - The vortex cores are denoted by '+' in both fig- 
ures, shown in the x-y plane. The maximum displacements 
of vortex cores are shown by the arrow. The maximal region 
around the core in (a) is elongated by the current. The simi- 
lar horizontal broadening around the core in (b) is caused by 
the vortex motion. Energy is transported; maxima in (a) and 
(b) do not coincide. 

In Fig.HJb), one can see that the system dissipates en- 
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ergy via vortex cores. From a microscopic point of view, 
Cooper pairs break into quasiparticles inside the core; 
these couple to the crystal lattice through phonons and 
impurities to transfer heat. The interaction between vor- 
tices and excitation of vortex cores manifests as friction 7 . 

The power loss of the system averaged over time and 
space is (P) = (Q). 



(P) 



e 

Jo 



d 2 



-e 2 



4oj 



(39) 

where I n is a Bessel function of the first kind. In Fig. [S] 



■ ujI\ 



Auj I 



0.01 



0.004 




FIG. 5: Power loss of supercurrent {J ■ e) (upper panel) and 
expansion parameter e (lower panel) as a function of frequency 



is shown the power loss and also e as a function of fre- 
quency, e is proportional the density of Cooper pairs, 
and can be thought of as an indication of how robust is 
the superconductivity. 

As frequency increases, while e tends to an asymp- 
totic value, (P) achieves a maximum and then decreases; 
this maximum is due to fluctuations of order parameter 
caused by the input signal. In a fully dissipative system 
as considered here, the maximum of each curve is not a 
resonance phenomenon but is instead caused by fluctua- 
tion of the order parameter resulting from the influence 
of the applied field. 

A parallel may be drawn between what we have ob- 
served in this section and the suppression of the super- 
conductivity by macroscopic thermal fluctuations com- 
monly observed in high-temperature superconductors. In 
our case, the vortices in a high-T c superconductor un- 
dergo oscillation due to the driving force of the external 
field. We may think of this as being analogous to the 
fluctuations of vortices due to thermal effects alone in a 
low-T c superconductor. Although the method of excita- 
tion is different, the external electromagnetic perturba- 
tion in the present case essentially plays the same role as 
the thermal fluctuations in low-T c situation. 

Finally, we point out that e seems to be an appro- 
priate parameter for determining the amount of power 



loss. Generically, it seems that for points deeped inside 
the superconducting region, that is at large e compared 
with its saturation value at high v, the power loss due 
to the dissipative effects of the vortex matter becomes 
suppressed. We suggest the possibility that this effect, 
which is naively intuitive, is in fact physical and more 
widely applicable than merely the present model. 



D. Generation of higher harmonics 

The practical application of superconducting materials 
is dependent on how well one can control the inherent 
nonlinear behaviour. In this section we will focus on the 
generation of higher harmonics in the mixed state, in 
response to a single-frequency input signal. 

The periodic transport current {J) r is an odd function 
of input signal, and it turns out that the response motion 
also contains only odd harmonics. From (|36|) we can 
calculate the Fourier expansion for transport current. 



i(2n+l)wTl 



n=0 



where the Fourier coefficient a( 2n+1 ') is 



,(2n+l) _ 



ee 2 



+ 1 



Auj 



-In 



AuZJ 



(40) 



-i(2n+l)9 



(41) 

We see the response goes beyond simple ohmic behaviour 
and the coefficients are proportional to e. Experimen- 
tally, one way of measuring these coefficients is a lock-in 
technique^ which is adept at extracting a signal with a 
known wave from even an extremely noisy environment. 

To make contact with more standard parameters and 
satisfy our intuition, we expand the first two harmonics 
in terms of v. The fundamental harmonic, expanded 
in powers of v 2 is 



r(l) 



ah 



0a(1 



iuo) 
ah 



(42) 



4/3a(1 
0{v A ) 



ioj) 1 



1 - Pb/Pa 



- or 



1 

ah 



i 

2uj 



where ah = (1 — t — b)/2b. The first term is the ohmic 
conductivity denoted as a\p , and is reminiscent of Drude 
conductivity for free charged particles. This is not an 
unexpected parallel, since the Cooper pairs in a super- 
conducting system can be imagined to behave like a free- 
particle gas. Taking this viewpoint, in the small-signal 
limit, the ratio Im cr/Re a — to = CIt s gives the relax- 
ation time of the charged particles. Subsequent higher- 
order corrections all contain uj in such a way that their 
contributions are suppressed at large uj. The coefficient 
of the n = 1 harmonic expanded in powers of v 2 is 



r(3) 



ah 



v 2 /lu 



8/3 A uj{Z - uj 2 ) ~ i{iuo 2 - I) 



0{v 4 



(43) 
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which decreases quickly with increasing uj. 

In Fig. [6l we show the generation of higher harmonics 
for three different states in the dynamical phase diagram. 
For each harmonic labeled by n, |cr' 2n+1 '| as a function 
of v has the same onset as e. We can see that |cr' 2n+1 '| 
reaches a maximum and then starts to decay while e sat- 
urates. The coefficients of harmonics with n > decay 
to zero in the high v limit, where the state is well inside 
the superconducting region. 

We pointed out in section lTlI Cl and reaffirm here that e 
plays a significant role in determining the extent of non- 
linearity in the system. In turn, the parameter which 
controls this is uj. When uj is large, e is brought closer 
to its saturation value eoo, causing the higher harmon- 
ics to be suppressed, and also lessening distortion of the 
vortex lattice. Finally, for a given harmonic, |er( 2 ™ +1 )| is 
generally smaller when is smaller; this can be seen by 
comparing (a) and (c) of Fig. [6l One might point out 
that the nonlinear behaviour is decreased at, for exam- 
ple, large uj. Nevertheless, we view the parameter e/e,^ 
as more intrinsic to the system, rather than simply char- 
acterising the input signal. 

A limited parallel can be drawn between the effect of 
thermal noise in high-T c superconducting systems, and 
the effect of the electromagnetic perturbation in our 
present case. It seems that in either case the fluctua- 
tion influence can be reduced by moving the state deeper 
inside the superconducting region. 



E. Flux-flow Hall effect 

In contrast to the fully dissipative system we have con- 
sidered, in this section we will discuss an effect caused 
by the non-dissipative component, namely the Hall ef- 
fect. In a clean system, vortices move without dissi- 
pation; a transverse electric field with respect to cur- 
rent appears. The non-dissipative part is subject to a 
Gross-Pitaevskii description, using a type of nonlinear 
Schrodinger equation 7 , with a non-dissipative part to the 
relaxation 7 from ([2|). The fully dissipative operator L 
in our previous discussion can be generalised by using a 
complex relaxation coefficient r = 1 + i£. We thus define 



C = rD T - l -{Dl + dl + d 2 z 



(44) 



The ratio ( = Im 7/Re 7 is typically on the order of 
I0~ 3 for a conventional superconductor, and 10 -2 for a 
high-T c superconductor—. 

The Hall Effect is small here. In normal metals, the 
non-dissipative part gives the cyclotron frequency. If r e 
is the relaxation time of a free electron in a dirty metal, 
then for typical values of uj c t c <C 1 the Hall effect be- 
comes negligible. Because the supply of conducting elec- 
trons is limited, the transverse component increases at 
the expense of the longitudinal component as the mean 
free path of excitations grows. It is equivalent to an in- 



crease in the imaginary part of the relaxation constant 
at the expense of the real part. 

The Eigenvalues and Eigenfunctions of C can be ob- 
tained easily by replacing the v in previous results with 
rv, uj with ruj and r with r/r. The transport current 
along the x-direction is no longer zero in the presence 
of the non-dissipative component; it is propotional to £. 
The frequency-dependent Hall conductivity can be ob- 
tained from the first-order expension in v, 



r Mi) 



c 



p A (1 - H 2 - w 2 C 2 



(45) 



while the Hall contribution in the y direction is expected 
to be negligible, as it is of the order of £ 2 . 

In principle, the crossover between non-dissipative sys- 
tems and dissipative systems can be tuned using the ratio 
£. In a non-dissipative system, which is the clean limit, 
the Hall effect is important and taking account of the 
imaginary part of TDGL is necessary. On the contrary, 
in a strong dissipative system where excitations are in 
thermal equilibrium via scattering, the TDGL equation 
gives satisfactory agreement. 



IV. EXPERIMENTAL COMPARISON 

Far-Infrared spectroscopy can be performed using 
monochromatic radiation which is pulsed at a high rate, 
known as Fast Far-Infrared Spectroscopy. This technique 
sports the advantage of avoiding overheating in the sys- 
tem, making it a very effective tool in observing the dy- 
namical response of vortices. In particular, one can study 
the imaginary part of conductivity contributed mainly 
from superconducting component. 

In Fig. [7] is shown a comparison with an NbN experi- 
ment measuring the imaginary part of conductivity. The 
sample has the gap energy 2A = 5.3 meV. The resulting 
value of 2A/T C is larger than the value expected from 
BCS theory 2 ^. We consider frequency-dependent con- 
ductivity in the case of linearly polarised incident light 
with a uniform magnetic field along the z axis. 

The theoretical conductivity contains both a supercon- 
ducting and a normal contribution. The total conductiv- 
ity is obtained from the total current as in where 
the normal-part conductivity in the condenstate is the 
conductivity appearing in the Drude model. 

According to our previous discussion, the nonlinear ef- 
fect of the input signal on NbN is unimportant in the 
THz region, which corresponds to v ~ 17. An approxi- 
mation where the flux-flow conductivity includes only the 
term ctq 1 ' from (l42j) . and Hall coefficient <Tq from (|45|) 
is shown in Fig. [7] and the agreement with experiment is 
good. 

The naive way in which we have treated the normal- 
part contribution is essentially inapplicable to the real- 
part conductivity. This is because the real-part conduc- 
tivity contains information about interactions with the 
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quasi-particles inside the core, making further consider- 



ation necessary^. 
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FIG. 7: Experimental comparison of imaginary part of con- 
ductivity at high frequency: The NbN experimental data are 
from Fig. 6(b) of H2]. Material parameters = 15.3 K and 
H" 2 = 14.1 T calculated using Fig. 3(b) of [gg]. Normal-state 
conductivity is a n = 20 • 10 4 (^if2 cm) -1 and relaxation time 
of an electron r e = 5 fs are taken from [2^ |. The theoretical 
curve has one fitting parameter n = 44.5. 



V. CONCLUSION 

The time-dependent Ginzburg-Landau equation has 
been solved analytically to study the dynamical response 
of the free vortex lattice. Based on the bifurcation 
method, which involves an expansion in the distance to 
the phase transition boundary, we obtained a perturba- 
tive solution to all orders. We studied the response of 
the vortex lattice in the flux-flow region just below the 
phase transition, at first order in this expansion. We 
have seen that there are certain parameters which can 
be tuned using the applied field and temperature, pro- 
viding a feasible superconducting system where one can 
study precise control of nonlinear phenomena in vortex 
matter. 

Under a perturbation by electromagnetic waves, the 
steady-state solution shows that there is a diamag- 
netic current circulating the vortex core, and a trans- 
port current parallel to the external electric field with 
a frequency-dependent phase shift and amplitude. Vor- 
tices move perpendicularly to the transport current and 
coherently. 

Using a technique of maximising the heat-generation 
rate, we showed that the preferred structure based on 
energy dissipation is a hexagonal lattice, with a certain 
level of distortion appearing as the signal is increased or 
the frequency is lowered. 

Energy flowing into the system via the applied field 
is dissipated through the vortex cores. We showed that 
the superconducting part may be thought of as having 
inductance in space and time. 
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We have written transport current beyond a simple 
linear expression. A comparison between different har- 
monics of three different states in our four-dimensional 
parameter space indicated that the nonlinearity becomes 
unimportant at high frequency and small amplitude, and 
the influence of the input signal is decreased when the 
system moves deeper inside the superconducting region, 
away from the phase-transition boundary. 

To observe the configuration of moving vortices, tech- 
niques such as muon-spin rotation2iI2, SANs^, STM 24 
and others 2 ^ seem to be promising options. To provide 
the kind of input signal considered here, methods such 
as short-pulse FIR Spectroscopy as used in [22| might 
be applied. The coefficient we defined in (|4Tj) corre- 
sponds to conductivity. We have also seen that a sim- 
ple parametrisation by complex quantities like conduc- 
tivity and surface impedence is insufficient to capture the 
detailed behaviour of the system; in performing experi- 
ments, it should be kept in mind that the nonlinearity 
can be measured in terms of more appropriate variables 
as we have shown. 

We have viewed the forcing of the system by the ap- 
plied field to be somewhat analogous to thermal fluctu- 
ations, in the sense that they both result in vibration of 
the vortex lattice. Hence, the influence of the electro- 
magnetic fluctuation is stronger at the nucleation region 
of superconductivity than deep inside the superconduct- 
ing phase. Besides, since at high frequency the motion 
of vortices is limited, the influence from electric field is 
suppressed, as is the Hall effect. 
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Appendix A: TDGL parameters 

We follow [3] to estimate the coefficient 7 which char- 
acterises the relaxation process of the order parameter. 
7 is the inverse of the diffusion coefficient for electrons in 
the normal state. For a strongly-scattering system, as in 
the dirty limiti, the ratio between the relaxation times 
of order parameter 



2m*JaT°(l-t)\ 
and the vector potential (or current) 

P(7 n 2m* h 



(Al) 



TA 



14C(3) 



(A3) 



By definition of the thermal critical field H c and H C 2 = 
as in Q , we know the ratio of the two parameters 



is 



(aT°) 2 (ff° 2 ) s 



(A4) 



The coherence length at zero temperature can be writ- 
ten in terms of H® 2 as £ 2 = $o/27riJ" 2 and in terms of 
effective mass as £ 2 = fi 2 /2m*aT c °. As a result, 



2w 5 K 2 (T r , 



(A5) 



With 7, we can retrieve the experimental quantities from 
the rescaled ones used in calculation. Using the 0- 
subscriptcd original variables, we write the electric field 
)e, and the frequency fi = (2/7£ 2 )ia 



E = (2h/e*£ 3 j) 
The current density Jo 



In the 

case of linear response, we have Jo = foE where 00 = 

(c 2 7/47TK 2 )<7. 



Appendix B: Solving the linearised TDGL equation 

We consider the linearised time-dependent Ginzburg- 
Landau equation ([9]) , which has been written in our cho- 
sen gauge. We wish to find the set of Eigenfunctions 
of L corresponding to the lowest Eigenvalue. Based on 
knowledge of the solution in the static case, we solve 
the now time-dependent problem by making the follow- 
ing Ansatz 2 ?' 2 ?. The electric field along the y-direction 
breaks rotational symmetry in the x-y plane, so we write 

f(x,y,z 7 T) = e lk * z e lk * x e 92< - T)y2+9l{T)y+9a{T) . (Bl) 

After substitution of / for $ in © , comparison of co- 
efficients of powers of y gives the following differential 
equations in r. 
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g\ + iv coswt — k x — 2g 2 9\ = 
90 + ^(k% + k 2 z ) - ig 2 - g 2 = 



0. 
0. 

£. 



The solutions are 
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-— tanh[c2 + t\. 



(B2) 
(B3) 
(B4) 

(B5) 



For a steady-state solution, r — ¥ 00, we have gi 
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As with gx, we have here C2 = 0. 



.g = -ik x — sin(wr 



— sin 2 (ljt 
8uj V 



c (B7) 



where Co is a normalisation constant. 
The resulting Eigenvalue is 



2 



4(1 + w 2 ) 



(B8) 



Now, although in a more realistic treatment of the sys- 
tem, one may introduce some boundary condition re- 
stricting k z to a certain set of values, here we simply se- 
lect the smallest Eigenvalue Eq available to us by setting 
k z to zero. Thus equipped with the set of Eigenfunctions 
corresponding to our lowest Eigenvalue, we deem them 
to be the first elements of our basis, labeled by n = 
and k x . These Eigenfunctions of L are 



<Po^(x,y,z,r) = e ik ^- vsinu ^u kx (y,r) (B9) 



with 



Uk{y,r) = c(r)e 



— \y~ k x -\-iv cos(cjt — #)J 



and 



c(r) 



^ [sin 2 0+cos2(u-r-fl) + ^ sin2(uT-fi)] 



(BIO) 



(Bll) 



In the same way, the corresponding Eigenfunctions of 
U can be obtained. 



<p , kx (x,y,z,r) =e ik * ( - x - vsi ™ r ^u kx (y,T) (B12) 



with 



coslcjr 



(B13) 



where 



c(r) 



_ p — ^ [ sin2 S+cos2(wT+fl)- i sin2(uT+0))] 



(B14) 



according to our normalisation condition 
est Eigenvalue of ]J is £q = £o- 



The low- 



In the ui — > limit, the system reduces to the case of 
constant electric field. The Eigenfunctions and Eigenval- 
ues are then consistent with those obtained by Hu and 
Thomson^. In the limit of zero electric field, the Eigen- 
functions and Eigenvalues reduce to those of the lowest 
Landau level static-state solution^. This is also the same 
in the u> — > oo limit. 
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